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ABSTRACT 



5-H , We investigate the formation of spherical cosmological structures following both dark 

^5^1 matter and gas components. We focus on the dynamical aspect of the collapse assuming an 

' adiabatic, 7 — 5/3, fully ionized primordial plasma. We use for that purpose a fully Lagrangian 

. hydrodynamical code designed to describe highly compressible flows in spherical geometry. 

We investigate also a "fluid approach" to describe the mean physical quantities of the dark 
^ ■ matter flow. We test its validity for a wide range of initial density contrast. We show that 

an homogeneous isentropic core forms in the gas distribution, surrounded by a self-similar 
I hydrostatic halo, with much higher entropy generated by shock dissipation. We derive analytical 

■ expressions for the size, density and temperature of the core, as well as for the surrounding halo. 

, We show that, unless very efficient heating processes occur in the intergalactic medium, we are 

unable to reproduce within adiabatic models the typical core sizes in X-ray clusters. We also 
show that, for dynamical reasons only, the gas distribution is naturally antibiased relative to the 
total mass distribution, without invoking any reheating processes. This could explain why the 
gas fraction increases with radius in very large X-ray clusters. As a preparation for the next 
study devoted to the thermodynamical aspect of the collapse, we investigate the initial entropy 
' level required to solve the core problem in X-ray clusters. 



Subject headings: Cosmology: theory - dark matter - galaxy formation - hydrodynamics 
methods: numerical 



1. INTRODUCTION 

This paper is the first in a series devoted to the hydrodynamical study of the formation of structures 
in cosmology. We focus here on a precise description of the dynamics of dark matter and gas collapse. 
Non adiabatic processes, namely temperature decoupling between the various plasma components, heat 
conduction, non equilibrium ionization balance, cooling, heating by star formation are examined in a 



following paper (Chicze, Teyssier & Alimi 1997) 



Gas and dark matter distributions observed in X-ray clusters provides strong constraints on structure 
formation theory. Although the gas density profile shows in the most massive clusters a typical core-halo 
structure (Hughes 1989; Briel, Henry & Bohringer 1992), the total mass distribution deduced from Icnsing 
analysis seems to be strongly peaked towards the center, with no clear evidence of a central core. EINSTEIN 
and ROSAT observations of X-ray clusters show gas core radii with sizes ranging from 100 to 300 h^^ kpc. 
Recent studies of gravitational lensing in clusters suggest that the total mass distribution can be singular 
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in the center (Wu & Hammer 1993 and reference therein) or at least with very small core radii {rcore — 
20-30 h-i kpc) (Tyson, Valdez & Wenk 1990; Mellier et al. 1994; Flores & Primack 1996). Moreover, 
gas and dark matter distribution differ also in the outer parts of galaxy clusters. David, Jones & Forman 
(1995) (hereafter DJF) show that there is a strong antihias of the X-ray gas in the clusters with respect to 
the underlying total mass distributions. The gas mass fraction increases with radius, up to what should be 
the primordial gas fraction VLb/^q — 0.2 — 0.3 (White et al. 1993). The origin of such an antibias is an 
open question. A model of cluster formation, which takes into account the dynamics of both dark and gas 
component, should be able to reproduces these features. 

The most rigorous approach so far of the coupled gas and dark matter dynamics was performed by 
Bertschinger (1985) in the case of the self-similar secondary infall. Its analysis showed that the final gas 
distribution traces the dark matter distribution. For more general initial conditions, it is generally admitted 
that such a conclusion subsists as long as thermal processes are not included. More recent studies, based 
essentially on 3D numerical simulations of structure formation, show however that gas and dark matter 
have different distributions due to the hierarchical nature of gravitational clustering. Antibiasing is due to 
a systematic energy transfer from dark matter to gas during mergers of protoclusters (Navarro & White 
1993; Pearce, Thomas & Couchman 1994). We propose to address again this question, by considering first 
the role of the purely adiabatic dynamics (this paper), and second the role of various thermal processes 



(Chieze, Teyssier & Alimi 1997) with a high resolution hydrodynamical treatment. 



The hierarchical nature of gravitational collapse requires, in principle, high resolution, 3D simulations. 
But it is difficult to reconcile physical completeness with topological realism. Plain hydrodynamics 
simulations generally suffer from limited resolution and reduced physical inputs, while one-dimensional 
calculations offer a practically unlimited resolution, with however a highly constrained and unrealistic 
geometry. They are also useful tools to work out and reveal the potential effects of various physical 
processes. Moreover, accurate numerical methods can be more easily refined in ID models, especially in a 
high resolution Lagrangian framework, and then provides a valuable complement to more elaborated 3D 
simulations of the formation of large-scale structures. 

The object of this paper is to study the gas and dark matter dynamics, with a one dimensional, high 
resolution, Lagrangian model, in a way similar to Thoul & Weinberg (1995) in spherical geometry or to 
Shapiro & Struck-Marcell (1984) in slab geometry. In the spherical case presented in this paper, we focus 
on the differences between gas and dark matter distributions due to dynamical processes only. In order 
to understand clearly the dynamical features of both components, it appears more convenient to consider 
separately the pure gas and the pure dark matter collapse, before investigating the coupled dynamics. 

When the calculation of the pure gas collapse of a perturbed spherical system is worked out along 
these lines, it appears that the collapse of the inner regions actually proceeds to very high densities. This 
is due to the fact that - given any reasonable density perturbation - the inner subsonic flow is not far from 
homologous, so that shock heating has here little effects. More precisely, it exists a central core which 
evolves strictly adiabatically up to hydrostatic equilibrium. This core is bounded by the sonic radius (where 
the gas velocity is equal to the sound velocity), from which the accretion shock originates and propagates 
outward. The final central density depends on the core initial entropy. We develop in Section 2.2 an 
analytical model which describes precisely the core structure and the surrounding halo. Section 2.3 presents 
results of numerical simulations which confirm both the numerical model developed in section 2.1 and the 
previous analytical model of a pure gas collapse. 

The collapse of collisionless dark matter is supposed to provide the potential well in which gas is falling. 
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and thus deserves special attention. The very high dynamical range required for the treatment of the gas 
component claims for an equivalent dynamical range in the treatment of the collisionless component. A 
classical method, adapted to spherical geometry, consists in calculating the trajectories of concentric shells 
of collisionless material. However, it appears rapidly that the numerical accuracy which can be achieved 
with this method cannot, for any reasonable computing time, meet the very stringent dynamical range 
and resolution requirements imposed by the gas evolution. First, the dark matter distribution in each 
Lagrangian gas element must be well sampled by a sufhcient number of dark matter shells, which requires 
a very large number of shells, with no guarantee of a smooth dark matter distribution. Second, the mass 
distribution must be calculated over more than 10 decades in radius. In view of these difhculties, we 
propose in Section 3 a "fluid model" for the dark matter component, which very accurately accounts for 
its gravitational behavior. In this way, the coUisional and collisionless components are treated on an equal 
numerical footing. 

The resulting coupled model contains all the dynamical aspects of the spherical collapse of both gas 
and dark matter. Some salient features of this model are discussed in Section 4. We show that for purely 
dynamical reasons the gas and dark matter distributions differ significantly, leading to an antibiased gas 
distribution relative to the total gravitational mass. We investigate the influence of the baryons density 
parameter on the gas distribution and we compare our results with the observational data of DJF, finding 
excellent agreements between our model and the most massive X-ray clusters. The difficulty to form 
realistic cores, even in adiabatic conditions, are finally discussed, in order to evaluate the initial entropy 
level required to explain the observed sizes of X-ray clusters core radii. It turned out that photo-ionization 
processes are not sufficient, and that galactic winds are good candidates to transfer enough energy to the 
intergalactic medium and solve the problem. The more complex and richer thermodynamical aspects of 



cosmic structure formation are discussed in a forthcoming paper (Chieze, Teyssier & Alimi 1997). 



2. GAS DYNAMICS 

In this section, we analyze the basic properties of the dissipative heating of the gas during gravitational 
collapse. The adiabatic dynamics does not depend on the ionization degree of matter. The gas is here 
assumed to be fully ionized, composed of 76% hydrogen and 24% helium by mass. We consider here the 
gravitational collapse of a purely gaseous, gravitationally bound, sphere. Emphasis is put on the treatment 
of shock waves, especially on the conditions of their appearance, both in time and space. Our high 
resolution Lagrangian treatment puts forward the formation of an adiabatic - isentropic - high density core, 
in which shocks never occurred. The importance of this core is that it determines the rest of the structure, 
essentially heated up by the dissipation of supersonic infall motions, in the relaxation layer of a moving 
accretion shock. 

Mass motion kinetic energy is dissipated into heat across this shock by the action of viscous forces. 
The dissipative shock region has a characteristic width of a few particle mean free paths, generally well 
beyond the resolution of an hydrodynamical simulation. In most treatments, shock transitions are handled 
with a so-called "artificial viscosity" by which the entropy production rate is largely increased above its 
actual value, so that the shock width is spread over few cells. The success of the operation relies on the 
fact that Hugoniot-Rankine relations, which connect upstream and downstream gas properties (density, 
momentum and energy) are conservation relations, independent of the dissipative processes which drive the 
transition. However, it is essential to respect the correct dependency of the viscous force on the structure 
of the velocity field, in order to determine safely where and when in the flow "artificial" viscosity is 
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required. The viscous forces are described by the stress tensor, which, in spherical geometry, depends only 
on the combination (§7 — -p)- As required, the stress vanishes for locally isotropic flows. This is of great 
importance for cosmological flows, which spend a long time close to uniform expansion or contraction. 

We present briefly the main features of the numerical code we employed for this study. It is an evolution 



of a code which has been developed for various astrophysical purposes ( Blottiau, Chieze fc Bouquet 1988 
de Boisanger, Chieze k Meltz 199^ ; [Rothcnflug et al. 19941 ). 



In order to get very high resolution, the numerical scheme is Lagrangian. In all calculations presented 
here, the number of Lagrangian shells is 250, both for gas and dark matter. Since the resolution must 
be concentrated in the core, shell masses are distributed according to a geometrical progression, with a 
dynamical range of 10^". For example, in the simulation of a 10^^ Mq cluster, the inner shell has a mass 
of 10* M0 only, while the outer shells have a mass of lO^'' Mq. It turns out that increasing the number of 
shells is unnecessary. The effective spatial resolution in the centre (in our case 0.1 pc) depends on the inner 
shell masses, not on the total number of shells, which just governs the resolution in the outer regions. 

Numerical stability is guaranteed by a fully implicit method, which allows to retain the time step which 
limits the relative variations of the dependent variables to some small arbitrary value, generally 1% in the 
simulations presented here. The independent variables are the time t, and the cumulative mass m from 
the center. The dynamical dependent variables are the radius r(m, t) of the sphere containing the mass 
m and the corresponding velocity u{m,t). Dealing with perfect gases of adiabatic exponent 7 = 5/3, the 
thermodynamical variable we adopt is the function 

s{m,t) = In (^EV^^ (1) 

where E is the internal energy of a fluid element of volume V. The evolution of this thermodynamical 
function, reminiscent of the entropy, is given by 

- = ^ (2) 
dt E ^ ^ 

where Q is the variation rate of the internal energy due to dissipation or other thermal processes, excluding 
the work done by the thermal pressure. A choice of this nature is dictated by energy conservation 
considerations : it is clear that isentropic flows (Q = 0) numerically remain strictly isentropic, since we 
avoid the integration of the work done by pressure forces over the whole history of the system, which 
invariably introduces numerical inaccuracy. The internal boundary conditions are naturally r (0, t) — and 
u (0, t) = (no artificial core radius). Since the perturbation we study are always compensated, the velocity 
at the outer boundary is taken to be the unperturbed Hubble flow. Our study is devoted to adiabatic flows. 
Therefore, the dissipation term Q is due to shock heating only. We now describe the form one should use 
for Q in spherical geometry. 



2.1. Dissipation in Shock Waves 

We choose to handle shocks by the pseudo-viscosity method (PV method) , introduced for planar shocks 
by Von Neumann & Richtmyer (1950). Though it is often considered as a numerical trick, this method has 
in fact firm physical grounds. In the PV method, the particle cross-section is reduced so that their mean 
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free-path matches the cell size of the simulation. Consequently, the PV method consists in solving correctly 
the Navier-Stokes (NS) equations with the available resolution, adapted in fact to a highly viscous flow. 
Note that the Hugoniot-Rankine relations, which connect the downstream to the upstream regions of an 
adiabatic shock, are independent of the detailed microphysics, since they just express conservation laws. 
This is the reason why the shock adiabats, calculated by the PV method, are strictly identical to those 
obtained by a NS calculation (see Chize, Teyssier & Alimi 1997 for further comments). 

In a fluid element, two distinct groups of time scales can be derived, related respectively to the 
microscopic properties of the plasma and to the macroscopic properties of the flow. In the first group, we 
find in particular the mean collision time TcoU between particles, and the relaxation time Trei of the velocity 
distribution due to collisions. In the second group, two time scales characterize the rate at which the shape 
of the fluid element changes, and the rate at which its volume varies. The ratios of the time scales of 
different groups, T(-]^-)/t(2), measure the departure from thermodynamical equilibrium, which results in the 
apparition of a viscous stress of magnitude 11 ~ p T{i)/'''{2)j where p is the thermodynamical pressure. 

Departures from thermodynamical equilibrium resulting from the finite relaxation time of internal 
degrees of freedom, give rise to the so-called bulk viscosity. The corresponding scalar correction to the 
thermodynamical pressure is proportional to V • m. Since the relaxation time of the energy levels of atoms 
and ions (~ 10~^s) is negligible compared to any relevant dynamical time scale pertinent to the problem, 
the bulk viscosity is absent. 

The usual viscous stress tensor is related to the local rate of strain, that is, the rate at which the shape 
of a fiuid element is altered. This means that locally isotropic motions can not generate any viscous stress. 



In Cartesian coordinates, the viscosity stress tensor, aij, and the dissipation rate have the form ( Landau 
Lifshitz 198^) 



, dui duj 2 ^ „ \ , A dui 

where is the viscosity coefficient, proportional to the particles mean free-path. Since we are interested 
here in spherically symmetric flows, we recast these expressions into spherical coordinates 

„ fdu 1 \ , 1 
= 2r/ I — - -V -M I and a,, = ct^^ = -- ct,,^ (4) 

The equation of motion has the form 

du dP 9$ da 3 

P-n = — ^ + y- -^^r (6) 

at or or or r 

where du/dt is the material acceleration of the fluid. This formulation of viscosity is of common use in 
the study of highly non- homologous flows (e.g. Tsharnuter & Winkler 1979; Cioffi, McKee & Bertschinger 
1988). 

Compared to the above NS solution, the single manipulation in the PV method concerns the viscosity 
coefficient r], which is used in the simulation. Since the actual particle mean free-path I is replaced by the 
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local cell size, Ax, the viscosity coefficient is increased by the factor Ax /I. Note that the radial component 
of the stress (eq. [Q) can be either positive or negative, which represents respectively a "tension" or a 
"pressure" . Since shocks lead always to compression, viscous forces act as a pressure. Artificial viscosity 
should then operate only when the velocity field obeys the conditions 

3^ < V • tt < (7) 

or 

This completes the prescription that we use to handle shock waves and to avoid spurious viscous effects. 



2.2. Pure Gas Adiabatic Collapse 

Let us consider first an adiabatic gas collapse (Slo = f^s = 1). Note that our results, obtained here 
in an Einstein-de Sitter universe, are also valid for different models, provided that the appropriate change 
of the cosmic time are made. We show analytically and we confirm numerically that a high density core 
forms, where no shock dissipation has occurred, surrounded by a self-similar hydrostatic halo, where shock 
dissipation is very strong. 



2.2.1. Initial Conditions 

We start at an early epoch i^, into the linear regime. At that time, every spherical shell has an 
initial radius r^. We assume that the velocity field is unperturbed, relative to the homogeneous expansion 
Vi — H{ti)ri, and that the initial density contrast has the form 

/ „ sinoX , , 

P« = pJl + <5o^j (8) 

This density perturbation does not correspond to a pure growing mode. As long as the initial time is choosen 
early enough, the decaying mode vanishes quickly. The Lagrangian coordinate of each spherical shell is 
defined as g = ^s^i/Rs where tan4>s = $s {^s — 1.437r). The total mass of the system (M = 4/37rpii?|) 
is then unperturbed. The collapse epoch of the system is defined as the first shell-crossing time, which 
writes for these initial conditions 

_ 3^1+^0 . . 

A r3/2 
"0 

The two independent parameters which fix then the history of the system are the total mass M and the 
collapse epoch tc. The initial uniform temperature of the gas is given by the cosmic radiation background, 
which is tightly coupled to matter from the recombination epoch at z = 1000 up to z = 200. This leads for 
Zi < 200 (see e.g. Anninos & Norman 1994) 
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In Section 4, we allow for much higher temperature due to possible reionization processes by stars or 
quasars UV background radiation fields. 



2.2.2. Core Formation 



As long as the pressure gradient term remains small compared to the gravitational term (eq. j^]), the 
pressureless solution of spherical collapse can be used. This regime is valid until t ~ tc. At that time, it is 



possible to derive precollapse analytical solutions ( Moutarde et al. 1995 ) . With the notations we use in this 
paper, this leads to the density profile 



(11) 



This behavior is obtained for each initial density profile for which dpi/dr — and d'^pi/dr'^ ^ at r-i — 0. 
Note that such a power law is valid for (7 <C 1 and t ~tc. A useful result is the position of the Lagrangian 
fluid element q at t ~ t^. 



Rs 1 



(12) 



One therefore deduces the precollapse density profile p (x r During this pressureless infall, each fluid 

element is adiabatically heated to a temperature given by 



T{q,t)^T, 



p{U 



(13) 



which allows us to compute the local sound speed for each fluid element Cs{q,tc) oc T{q,tc)'^/'^ . As the 
velocity field of each infalling fluid element is 



u{q, tc 



Rs 



TT 1 
-g3 



4$5'5o 2tc 

we can deduce the Mach number of each fluid element defined as M{q, tc) = —ujcg 



(14) 



>((g,i.) = (y 



Rs 



27r 



4$s<Jo tcCs(tc 



5 

g3 



(15) 



The fluid element for which A^(g,tc) = 1 reaches its sonic radius at t ~ tc- Its coordinate is of primary 
importance in our discussion. Let us call the region enclosed by this particular shell the core of the forming 
cluster. For the fluid elements with q < qcorei the flow is always subsonic, and no shock dissipation occurs 
in this region. The flow leads to an hydrostatic configuration with nearly uniform density, pressure and 
temperature. 

Outside the core, the flow is supersonic, with higher and higher Mach number for increasing radii. 
When shocks are strong enough, we enter another regime of the flow leading to a self-similar structure we 
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identify as the halo. Its features are discussed in the next section. Between the core and the halo, there is 
an intermediate region where shocks are weak and heating is a mixture of adiabatic compression and shock 
dissipation. This ensures the transition between the adiabatic regime (core), and the strong shock regime 
(halo). Because the core is in hydrostatic equihbrium, its features remain remarkably stable up to z = 0. 

It is possible to derive analytically within this adiabatic model the core density, temperature and 
radius. Using equation (|l5|), we deduce the value of Qcore- We then obtain with equation ( [l^ ) the core 
radius of the cluster 



Tq is the average temperature of the unperturbed background at collapse epoch, adiabatically extrapolated 
until today. M is the cluster mass and Zc the collapse redshift. Note that the core radius is very sensitive 
to the temperature. The core entropy is indeed equal to the mean background entropy. Its size is directly 
related to the chosen value of the initial entropy. Because the average baryon density in a purely baryonic 
universe is fixed by the cosmological parameters, this means that the core radius is mainly determined by 
the value of the temperature just prior to the collapse. The final density inside the core can be computed 
using equation ( pT| ) 

_.,.,„w.(l±i)*(^)^(^)-^ ™ 

This core density depends strongly on the collapse epoch, because of the (1 + z)^ dependence of the 
mean background density. The core temperature is simply derived using equation (|l 

T_ . 8 X 10«K (l±^) " ^ f ) (18) 

These scaling laws agree quite well with the numerical resolution of the hydrodynamics equations, as 
we will see below. 



2.2.3. Halo Formation 

As soon as the core has been formed, Lagrangian fluid elements fall towards the center with supersonic 
velocity, leading to strong shock heating. Let us suppose that a given fluid element is shocked at an epoch 
given by its own collapse time, in the limit g <C 1, this gives 

ts{<l) = tc{l + ^q^^ (19) 
from which the Lagrangian coordinate of the shock front, at any time after the collapse, is obtained as 



(20) 
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Using equation (|l2|), we obtain tlic time dependence for the shock radius, which scales as cc {t — tcY^^. 
As soon as the shock wave is strong enough {Ai ^ 1), it is possible to estimate the postshock temperature 
T2 and density n2 using Hugoniot-Rankine relations. If we assume that the preshock temperature Ti 
is negligible (this is true by definition for a strong shock) and that the postshock velocity U2 vanishes 
(hydrostatic flow), this leads to the results 

kT2 = \l-i.mHu'^{qs^tc) and 77.2 = 4^^^^^-^ (21) 

After the shock front, because 1*2(9) = Oj density and temperature remain constant in time. This 
implies that the density profile after shock dissipation has the same slope as the precollapse density profile 
p oc r~^^/^ (in section 3, we will see that for dark matter the picture is somewhat different) and that 
the temperature profile scales as T oc r^/'', slowly increasing with radius. As the accretion shock escapes 
towards the outer part of the cluster, the small q approximation breaks down and the halo features are 
much more complex. 

The basic approximation we make here is the assumption of an hydrostatic postshock flow. When 
we computed the core properties, we have also neglected all pressure effects which can modify the 
trajectories. These points have to be tested with numerical calculations which allow us to solve precisely 
the hydrodynamics equations. 

2.3. Numerical Calculations 

To test all the formulae presented in the upper sections, and also to illustrate the importance of 
the numerical prescription used in this paper to treat shocks, we run a simulation with the same initial 
conditions as previously stated. The run parameters in this section are 

M = 10^^ Mq Zc^2 To = 0.014 K VLb^I (22) 

We start the simulation at Zi — 200 and reach z = with approximately 10^ time steps. We plot in Figure |l| 
density, temperature, velocity and entropy profiles obtained at different redshifts, using equation (^) for the 
viscous dissipation term. The density and temperature jumps at the shock front are just what is predicted 
by the Hugoniot-Rankine relations. The assumption of hydrostatic postshock flow appears to be valid with 
great accuracy, since the postshock velocity field is zero everywhere, and the density and temperature 
profiles remains constant in time. Note that no spurious oscillations appears in the post shock velocity 
field. The power laws for the density and the temperature profiles in the halo are in perfect agreement 
with the analytical derivation. By the virtue of the prescription for viscosity (eq. |^), the core remains 
strictly isentropic during the collapse (see fig. ||l|). Any other prescription would have led to a spurious 
core heating, making impossible any description of the most inner part of the cluster. 

More interesting is the good agreement between the computed values of the density and the 
temperature of the core and the predicted values. The computed core density and temperature are 
respectively 5 x 10''cm~'^ and 6 x lO^K, while the predicted ones were 7 x lO^cm"'^ and 8 x lO^K. The 
computed core radius is about 2 x 10"'^ kpc while the predicted core radius was 7 x 10~^ kpc. The 
difference is larger for the core radius, since pressure gradients, which have been neglected in our analytical 
approach, slow down infalling shells slightly before collapse. The important quantity here is the core 
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entropy, determined quite accurately by both the analytical and the numerical calculations. Note that the 
very small size of the core (~ 1 pc) has no influence on the large scale properties of the flow. 

When the entropy profile starts to rise, wc enter the halo regime, where strong shock dissipation takes 
place. When the shock front reaches outer regions of the cluster, the density profile get sharper, since the 
small q approximation is not valid anymore. The density profile slopes from a r"^^/'' law in the center to a 

law in the outer parts. 

3. A FLUID MODEL FOR DARK MATTER 

The coUisionless component of the universe is supposed to dominate the total gravitational potential 
at nearly all scales. In a high dynamical range simulation, the main difficulty is to get for dark matter a 
resolution comparable to that achieved for the gas component, ranging up to 10 orders of magnitude in 
density. 

We first approach the problem of dark matter dynamics by the traditional "spherical shells model" 
(Hoffman, Olson & Salpeter 1980; Peebles 1982; Hoffman, Salpeter & Wasserman 1983; Hausman, Olson 
& Roth 1983; Martel & Wasserman 1990; Thoul & Weinberg 1995). In this description, the coUisionless 
particles are distributed over A'' concentric shells of individual masses m and radius r. As long as the 
mass M included below a given shell is constant, its dynamics is determined by Kepler's laws for a central 
potential motion. The highest accuracy on kinematics and energy conservation (~ 10~^^ for a SUN station) 
is achieved by using this exact solution for each shell, up to shell crossing. At that time, new total energies 
can be calculated, and new exact trajectories are computed, until the next shell crossing time, and so on. 
This time-consuming procedure guarantees in principle a very high level of accuracy, since trajectories 
are analytically treated. However, in the limit of a continuous dark matter distribution with zero angular 
momentum, shell crossings occur continuously at the center of the spherical distribution, where the density 
is infinite. With a finite shell number, the first crossings (in radii) occur in a rapidly shrinking region, 
the size of which ultimately excludes any numerical solution. Moreover, a correct sampling of the shell 
distribution in the inner regions (in particular over the size of the adiabatic core of the gas) requires an 
excessive number of shells, resulting in a definitely prohibitive computing time. 

We have also examined various adaptations of the Particle Mesh (PM) technics, abandoning the quality 
of the semi-analytical solution. Again, the very small inner scales were not well-sampled, mainly because of 

the Eulerian nature of PM schemes. 

We propose in the following a fiuid treatment for dark matter, which reproduces properly some aspects 
of the self-similar solutions of dark matter collapse, derived in particular by Fillmore and Goldreich (1984) 
(hereafter FG) and Bertschinger (1985). This unusual choice, initially dictated by the high resolution 
requirements of Lagrangian simulations, provides a new insight on coUisionless dark matter dynamics. 
Moreover, it prevents for any two-body relaxation effect, inherent to any description of dark matter with a 
finite number of shells. 



3.1. Dgirk Matter in Spherical Geometry 

The coUisionless dark matter fluid is described by the velocity distribution function f{r,v,t), which 
evolves according to the Poisson-Vlasov equation. For a spherical system of particles moving along radial 
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orbits, it reduces to f{r,Vr,t), where Vr is the individual particle radial velocity. The hydrodynamical, bulk 
velocity of the fluid is u=< Vr >■ 

In a collapsing halo, three different regions can be distinguished: 

• The single stream, infall region, where the motion is purely Keplerian. The distribution function 
reduces to a 5 function, centered on the local flow velocity. 

• The relaxed core, where the dark matter fluid has been heated up by phase mixing. We show in 
a companion paper (Teyssier, Chieze & Alimi 1997) that, for the wide class of self-similar collapse, 

the adiabatic invariancc of the gravitational potential ensures in this region the evenness of the 
distribution function in velocity space. We take advantage of this property to derive the equation of 
state of the dark matter fluid. 

• The relaxation region, between the single stream and the relaxed regions, where the velocity 
distribution can only be described by solving the Poisson-Vlasov equation. We assume here that this 
dissipation region is of finite size, so that one can connect the flow properties on both sides, by the 
mass, momentum and energy conservation laws. The resulting "jump" conditions are the analogous of 
the Hugoniot-Rankine relations for classical fluids. Therefore, the preceding discussion of dissipation 
in fluids can be extended to dark matter. 



We first derive the equation of state of dark matter, in the relaxed region where our fluid approach 
is valid, from virial and from thermodynamical considerations (see also Teyssier, Chize & Alimi 1997 for 
further comments). 

We consider a large system of particles, with an equilibrium spherical distribution. Particles are moving 
only along radial paths, as required by the spherical dark matter model. Consider the N particles interior 
to the sphere (S) of radius L and mass M . The aforementioned evenness of the distribution fimction is 
equivalent to assume a quasi-stationary state, with no mass and energy flux through any surface element. 
An equal number of particles are flowing in and out the system, crossing (S) with opposite velocities. These 
two opposite flows result in a net momentum flux, pndS, through any surface element drawn on (S), with a 
direction opposite to the normal n to the surface element dS. It is thus without consequence to materialize 
the surface (S) by a reflecting sphere, on which the exterior fluid exerts a uniform, normal force —pn per 
unit surface. Defining the net force exerted on the surface (S) by the fluid exterior to the system by 



3.2. Dark Matter Fluid Thermodynamics 




(23) 



one can see immediately that the total kinetic energy of the system, E, obeys the virial equation 



2E ^ LuL 



(24) 



This is the equation of state of the system, which has been derived here from purely mechanical 
considerations, assuming a even particles distribution function. 
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The thermodynamical point of view leads to more specific results. The energy E and the entropy s are 
the state variables of the system. Both are extensive quantities relative to the total number of particles N. 
Since each particle has only one translational degree of freedom, the entropy, function of TV, E and of the 
external parameter L, has the form 

s {E, L) = Nk {\n ^ + i In ^) + (25) 
The first and second principle of thermodynamics lead to the expression 



dE = Tds - udL (26) 
where one has defined the intensive conjugated variables of S and L : 

T^C^) ^2 A and C.^-('^) =2^ (27) 
which is again the equation of state of the fluid. Clearly, for an isentropic evolution, one has 

dCj\ uj 

— ==3- 28 

where we have defined the linear mass density /i = dM / dL. It can be easily shown that (doj/dfi) g is just 
the square of the sound velocity. Specific heats are properly defined for transformations respectively at 
constant u) and constant L. From equations ( p^ ) and (|2^) one readily obtain the corresponding specific 
heats and Cl as 

= ^Nk Cl = ^Nk (29) 

With this peculiar understanding of the specific heats, this monodimensional gas can be considered 
as a 7'-^-' = Cq/Cl = 3 gas. The thermodynamical equilibrium of a (gravitational) force free system is 
characterized by the uniformity of T, uj which imply a uniform linear mass density ji. In other words, the 
equilibrium mass density p and particle velocity dispersion <j profiles are 



p oa L a = Const. (30) 

However, starting from an arbitrary distribution of particles in phase space, thermodynamical equilibrium 
can only be reached through gravitational interactions, since particles are supposed to be coUisionless. 



3.3. Dark Matter Fluid Equation of Motion and Dissipation 

In the dark matter spherical model, where particles have no angular momentum, no initial velocity 
dispersion, and suffer no collision, their orbits of zero eccentricity degenerate into purely radial trajectories. 
The motion of an individual particle is thus monodimensional, contrary to a classical gas for which collisions 
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populates the three translational degrees of freedom. Defining the hydrodynamical velocity u as the mean 
local velocity of particles, and the gravitational potential as ^(r, t), the equation governing adiabatic fluid 
motions is immediately derived from the definition of uj 



du 1 duj 9$ 
dt fi dr dr 

which express the conservation of the momentum of a spherical shell of unit mass. The position coordinate 
is now written as usual in spherical geometry r. The equation of motion (eq. is hyperbolic and thus 

leads to the formation of shocks. It cannot account for the phase mixing, which results, according to the 
Vlasov-Poisson approach, from shell crossings. As already stated, we assume that the relaxation region, 
identified as the first leading caustics in the solutions of the Vlasov equation, is of finite size. The upstream 
and downstream flow properties are connected by the continuity of the mass, momentum and energy flows, 
which write 



^2^2 = fJ-iUi (32) 

012 + t^2u\ = U)l+ filul (33) 
^U!2 + ^fl2ui] U2 = (^COi + ^^iujjui (34) 



where the subscripts refer to the upstream and downstream regions. Consequently, the "jump" conditions, 
which must be verified across the relaxation region, are 



^2 = 2^1 U2 = -ui UJ2 = -^fJ-iuj (35) 

Since they do not depend on the detailed relaxation processes, they are also the solutions, for our dark 
matter fluid, of the analogous of the Navier-Stokes equations (which however do not provide the exact 
solution inside the relaxation layer). Due to the symmetry of the system, the viscous forces acting on a 
fluid element is purely radial, so that the equivalent of the Navier-Stokes equation is 



du d , 9$ 

where u) is the viscous force defined on the surface of the sphere of radius r. We determine now the rate at 
which energy is dissipated in a fluid element of unit mass. The variation of its kinetic energy, can be 
written as 



dK d{ijju) d{ijju) 9$ _ du ^ du 

The first, second and third terms of the rhs represent respectively the rates at which the pressure, viscous 
and gravitational forces do work on the fluid element. This total amount of energy is not entirely restored 
as flow kinetic energy. Indeed, the two last terms represent the energy which goes into heat, the first of 
which is the compression work. Therefore, the expression 



represents the rate of increase of the internal energy resuhing from the dissipation of the kinetic energy by 
viscosity. Since dissipation can only result in an entropy increase, one has furthermore 




(39) 



The expression of Q must be definite and positive in dissipative regions, that is in shell crossing regions, 
near the leading caustic, and zero elsewhere. In this monodimensional (radial) fluid, there is no distinction 
between the rate of strain and the rate of increase of the length of a linear element of fluid, du/dr. For 
these reasons, the dissipative force can be written as 



which automatically verifies the irreversibility condition expressed by equation (^9|), provided that 77 > 0. 
This expression of the viscous force would be in principle valid for a coUisional monodimensional fluid, for 
which the particle velocity distribution is governed by collisions. This is not true for coUisionless material, 
where the particle velocity distribution remains completely determined by the initial conditions until the 
occurrence of the first shell crossing. Until that stage, the coUisionless fluid conserves the memory of the 
initial velocity distribution, contrary to a coUisional gas for which this memory is abolished after a few 
collision times. In our fluid description, it is reasonable to consider that shell crossing, and thus phase 
mixing, actually occurs when the criterion 



is fulfilled. We shall see in the following section that, according to this criterion, the position of the fluid 
relaxation layer (that is, the shock position) coincides quite accurately with the position of the leading 
caustic of coUisionless material. 



The role which is essentially assigned to dark matter is to provide the gravitational potential well in 
which gas eventually condenses. Therefore, two main requirements must be fulfilled by a fiuid treatment. 
At first, the over-all dark matter mass profile must be correctly described. Second, the time dependence of 
the growth of the relaxation region must be also correct. We check here these two points by comparing our 
fiuid model results to some known self-similar solutions of the problem of dark matter collapse. 

Similarity solutions have been found by FG for the problem of the collapse of cold, coUisionless matter 
in a perturbed Einstein-de Sitter universe. The scale-free initial mass perturbation they have considered is 
of the form 



du 



(40) 




3.4. Self Similarity Solutions for Dark Matter 
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where Mq is a scaling mass and < e < 1. The initial velocity field they have considered is unperturbed 
relative to the Hubble fiow. After shell crossing, the dark matter density profiles they obtained are well 
fitted by a power law p{r) oc r" at small radii, with, for spherical geometry, the following values for n 

2 9e 2 

n — —2 for e < - and n — for e > - (43) 

- 3 3e + 1 - 3 ^ ' 

We present in Figure (^ the overdensity profiles obtained at different redshifts for the cases e = 0.2, 
0.5, 0.8 and 1. The radius have been rescaled at each epoch to the self-similar radius defined as A = r/rta 
where rta is the current turn-around radius. Except naturally the caustics, the overdensity profiles we 
obtain are exactly those presented in FG. In our fluid model, dissipation occurs for dark matter through 
shock waves. The position of the shock front marks the first shell crossing, and therefore has to be close to 
the first caustic. We define the shock front as the point where the second derivative of the velocity field 
reaches its minimum. We plot the position of this shock front in Figure (||) as an arrow on the A axis. Note 
that both the power law and the position of the first caustic are in very good agreement between our fluid 
model and the coUisionless self-similar solutions. Note also that the dynamical range in radius we obtained 
is much larger than in any other dark matter simulations, allowing us to achieve very high resolutions. 

The case e — 1 corresponds to the secondary infall onto an initially overdense perturbation embedded 
in an Einstein-de Sitter universe. Bertschinger (1985) has examined this situation in details, both for 
coUisionless and coUisional fluids. Our result is quantitatively similar to the relevant coUisionless overdensity 
proflles presented in Table 4 of Bertschinger (1985). This author found that the position of the leading 
caustic for dark matter is A = 0.364. We recover this value with 2% accuracy within our fluid model. 
Bertschinger (1985) found that the shock front position for a, j — 5/3 gas was A — 0.339, very close to the 
dark matter leading caustic. He noticed however that 7 = 3, pertinent to a one dimensional gas, should have 
been preferred, although he obtained in that case A = 0.600. This apparent paradox has been elucidated 



here (section 3.2) 



3.5. Arbitrary Initial Conditions 

The similarity solutions for dark matter collapse are well recovered within our fluid model. It reproduces 
also the precollapse scale invariance obtained by Moutarde et al. (1995) for a wide range of initial density 
proflle. This precollapse scale invariance allows to built up naturally the conditions for self-similar collapse, 
for which our fluid model is fully valid. As a matter of fact, we plot in Figure (^) the evolution of the 
density profile for a pure dark matter collapse, using the same initial conditions as in section 2. Before the 
first caustic, the density profile follows the precollapse scale invariance p oc r~^^/^, and after shell crossings, 
the density is well fitted by the power law p oc r~^. It is important to notice that the mass distribution is 
not the same than for the pure gas collapse. Moreover, the postshock velocity vanishes exactly for the gas, 
leading to a constant central density, while for dark matter, the central density increases continuously with 
time, following the time scaling law given by FG for the case e = 4/9, which corresponds to a precollapse 
density profile with a slope n = 12/7. At late epoch, and for large radii, the self-similar regime is no longer 
valid, and a fully numerical treatment is necessary. 
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4. COUPLED DARK MATTER AND GAS EVOLUTION 

In section 2, we derived the evolution of a pure gas initial perturbation, while in section 3, we derived 
the evolution of a pure dark matter one, using for both case the same initial profile. We conclude that these 
two cases give different results for the mass distribution within the cluster (p oc r~^^/^ for gas and p oc 
for dark matter). In this section, the fluid we consider is a mixture of gas and dark matter particles. 
This leads to a new important parameter for the dynamics of the forming cluster, the gas fraction of the 
unperturbed background f^s (we still assume an Einstein-de Sitter universe with h = 0.5, but our results 
are valid for any cosmological parameters). 

4.1. Halo Structure 

For a dark matter dominated universe (with Qb ^ 1), it seems natural to obtain for the gas density 
profile very similar results than for dark matter. To test this hypothesis, we run three simulations with the 
same initial parameters 

M = 10^6 Mq = 2 To = 0.014 K (44) 
allowing the baryon density parameter to take three different values 

^3 = 0-01 f^B = 0.9 ris= 0.99 (45) 

We choose a low initial temperature, which, as we know from section 2, leads to a very small core 
{rcore — lO^'^ kpc). This is of course very far from the observed values in large X-ray clusters where 
rcore — 100 kpc. We are interested here in the halo structure, where a power law is observed for both dark 
matter and gas density profiles. The core problem will be addressed below. 

We plot in Figure (^) the gas and dark matter overdensity profiles obtained at z = for various Qb- 
For sake of comparison, we also plot the gas density profile for the fis = 1 case, and the dark matter 
density profile for the ils = case. Note that, even for low values of the baryon density parameter, the 
gas distribution does not follow the dark matter distribution. The gas density profile is steeper than for 
the pure gas case, but does not fit the dark matter profile. This is due to the presence of a core that 
prevents gas to fall deeper in the dark matter potential well, and to the shock relations followed by gas 
(three internal degrees of freedom) which differ from those followed by dark matter (one internal degree of 
freedom). Note also that the dark matter density profile is insensitive to the value of Og. The main effect 
of dark matter on the gas distribution is to compress the core and the halo after the collapse epoch. This 
leads to a strong heating of the core due to adiabatic compression only. The resulting temperature is much 
higher than for the pure gas case. The slope in the temperature profile changes also slightly, leading to a 
much more isothermal profile. For Qb — 1, the gas distribution recovers the same characteristics than for 
the pure gas case. 

One interesting quantity is the gas fraction, defined as fgas = Mgas{< r) / Mtot{'^ r), the ratio of the 
enclosed gas mass to the total gravitational mass. Because the slope of the dark matter density profile is 
steeper {pdark oc r~^) than for gas (pgas oc r~^^/^), the gas fraction increase with radius. In the core, fgas 
has a very low value, because dark matter dominates. In the halo, for a dark matter dominated Universe, 
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fgas increases continuously with radius as r^/^, with fgas — in the most outer regions. We plot in Figure 
(§) fgas versus the total mass as in DJF (see their fig. [5]). When we compare our profile with the observed 
ones in the most massive clusters, we find an excellent agreement. 

4.2. The Core Problem 

We saw in section 2 that, for very general initial conditions, a core forms with very high density, and 
we just saw in this section that if dark matter is present, the central density is even higher. The observed 
core density in massive X-ray clusters is not higher than 10~^ cm~'^, with core radii of the order of a few 
hundred kpc. This apparent discrepancy is due to the fact that in adiabatic models the initial entropy 
of the primordial gas is very low, leading to very small cores. We saw that within adiabatic models, it 
was impossible to generate entropy inside the core, even within a dark matter dominated scenario. When 
cooling is taken into account, the situation is worst: the resulting decrease in entropy naturally leads to 
the formation of an ultra high density central core, in about a free-fall time. It seems clear that this 
uncomfortable situation is unavoidable, if one restrict to dynamical considerations only. On the other 
hand, it is possible to imagine many entropy generating processes, such as supernovae driven winds, UV 
background radiation field, ambipolar diffusion, electronic conduction and so on. It is far beyond the scope 
of this paper to treat all this processes. To see an interesting review of this question, see Sarazin (1990). 
We simplify greatly the problem by considering a simple model, that we describe now. 

We consider an homogeneous heating process, such as a background UV radiation field generated 
by primordial quasars or stars. We model this process by turning on gas heating at say z = 5, and by 
imposing the following conditions. If T < 10^ in a cell, we reset the gas temperature by imposing T — 10^. 
On the contrary, if T > 10^, the gas temperature in the cell remains unchanged. This is a continuous 
entropy generating process, which has also physical support from the processes we just enumerated. In this 
framework, the background temperature is a free parameter, and we consider here 3 different values 

Tbg^lO^K Tbg^lO^K r.g^lO^K (46) 

The last value is much higher than the permitted temperature obtained by quasars spectra, but it gives 
us a limit for high entropy generating processes, such as supernovae driven winds. We take the following 
parameters for these three runs 

M = lO^*' M© = 2 VLb=0.1 (47) 
This method enables us to use the analytical formulae derived in section 2 for the pure gas collapse, using 

To = ng{i + z,)-^. 

We plot in Figure (||) the resulting structure we obtained at z = 0. Note that the core density is much 
lower than for the other cases. The core entropy has been greatly enhanced, leading, as expected, to a 
more extended core. Outside the core, the profiles remain unchanged by the heating process relative to the 
purely adiabatic case. Note also that we predict with our analytical model a decrease of the central density 
by a factor of 16 between the three runs, which is exactly what we recover in Figure (H). The dark matter 
overdensity profile is still insensitive to the gas distribution, as noticed in the last section. 

Note that even for the higher temperature run, the core size is still unrealistic. This is partly due to 
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the large mass and collapse redshift we use in our simulations. We believe however that in real clusters 
the initial entropy level should be as high as T^g ~ 10^ K. in order to obtain large enough core radii. An 
interesting study, that we postpone to a next paper, would be to deduce from observed X-ray clusters the 
required entropy level to explain their core sizes and densities. This could give important clues to seek for 
efficient heating processes in galaxy clusters. 

5. CONCLUSIONS 

The main goal of this paper is the study at very high resolution of the collapse of cosmological 
structures. This requirement implies necessarily a ID approach. Pancakes formation has been thoroughly 
studied over the past decade (Zcl'dovich 1970; Bond et al. 1984; Shapiro & Struck-MarccU 1985; Anninos 
& Norman 1994), while spherical collapse remained unexplored at an equivalent resolution (in our case, the 
central resolution is of the order of 0.1 pc). 

In the adiabatic case examined in this first paper, the main features of the virialized structures result 
from the dissipation of supersonic infall motions. One therefore should take care of the method used to 
handle shock waves, especially in spherical geometry. Within the artificial viscosity method, we derive the 
form of the viscous forces one should use to get a correct description of shocks, both in time and space. 
This general prescription can be easily adapted to 3D hydrodynamics. 

This provides a clear account of the formation of core-halo structures in the pure gas case. During 
the collapse, the flow is characterized by a subsonic inner region (the core), surrounded by a large envelop 
where strong shock dissipation occurs (the halo). The core and the halo are merely the remnants of the 
adiabatic and dissipative regions of the flow. The size of this isentropic core is therefore mainly determined 
by the initial entropy of the primordial gas. We obtain for a wide range of initial conditions very small cores 
with very high density, definitely ruled out by X-ray observations. 

These features are not qualitatively affected by the presence of the dark matter component, but 
quantitatively the situation is even worse. The main effect of dark matter is to compress further the gas 
core, but with subsonic velocities, so that, again, no entropy is generated. 

However, this conclusion relies on the treatment of the inner gravitational field dominated by dark 
matter. Due to the divergence of the potential at the center, methods using particles are excluded in 
practice by the current computers limitations. For that reason, we propose a fiuid approach describing the 
dynamics of dark matter. This method successfully reproduces the analytical solutions of the self similar 
spherical collapse. However, caustics are smoothed out within this model, but the leading shell crossing 
appears as a discontinuity in the dark matter flow at the correct location. This model provides moreover a 
new insight in dark matter spherical dynamics, treated here as a 7 = 3 mono dimensional fluid. 

We specifically studied initial conditions with a single spherical Fourier mode as initial density 
perturbation. The results we obtain are then valid for a wide range of arbitrary spherical initial conditions. 
We show analytically that in the pure gas case the halo has a density profile p oc r^^^/^, while in the 
pure dark matter case, it is close to p oc r~^. The numerical results are in perfect agreement with this 
scale invariance. The specific features of the coupled case are the followings. We show that, even for very 
low VLb, the gas never strictly follows the dark matter particles. This is due to their different numbers of 
translational degrees of freedom. As a consequence, the gas mass fraction is an increasing function of radius. 
Gas is more smoothly distributed than dark matter, with a density profile in-between those obtained for 
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pure gas and pure dark matter collapse. Adiabatic models are then an natural explanation of the observed 
antibias in X-ray clusters (DJF). In the most outer regions, the gas mass fraction reaches the unperturbed 

value ^g. 

These conclusions can be naturally affected by the actual geometry of a fully 3D gravitational dynamics, 
in which the collapse occurs in sheet-like structures, in filaments and also in quasi-spherical knots. This 
could lead to different density profiles, mainly in the inner part where rotational effects are expected to 
be important. Recent 3D studies show indeed that the mass profile in clusters may be smoother than the 
singular isothermal sphere in the central part (p oc r~^; see Navarro et al. 1995; Cole & Lacey 1996; Tormen, 
Bouchet & White 1996). However, even in 3D, the dense adiabatic cores are likely to merge subsonically, 
leaving unsolved the core problem. Anninos & Norman (1996) investigated the influence of the actual 
numerical resolution on a simulated cluster. They conclude that the gas central density does not converge, 
even for their highest resolution run. In fact, we show in this paper that the typical size of cores obtained 
for an adiabatic collapse should be about 10~^ kpc. Heating processes are therefore unavoidable to obtain 
measurable core radii. Antibiasing has been also investigated in 3D simulations (Kang et al. 1994; Martcl et 
al. 1994), where a strong segregation between gas and dark matter has been detected. However, Steinmetz 
& White (1996) have shown that this effect could be due to a spurious heating of the gas, due to the finite 
number of particles used in 3D simulations. We stress here that the antibias discussed in this paper, in a 
rather constrained geometry, has a firm physical origin: the different dimensionality of gas (7 = 5/3) and 
dark matter (7 = 3) phase spaces. This result was obtained using analytical or semi-analytical methods, 
free from any two-body relaxation or spurious viscous effects. 

We briefly investigate in this paper the initial entropy level which would lead to more realistic cores 
properties. We find that UV radiation fields, accounting only for a background temperature of 10^ K up to 
10^ K, are unable to solve the core problem. Consequently, it appears that the initial entropy level that 
could alleviate this short cut is very high, and must result from a background temperature as high as 10'' 
K. We therefore analyze in a companion paper the thermodynamical aspects of collapse, regulated by the 
balance between gas cooling and new energy sources like the feedback from star formation. 

We thank an anonymous referee for his careful reading of the manuscript and his constructive comments 
which have improved the quality of our work. 
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Fig. 1. — Density, temperature, velocity and entropy profiles obtained at z= 1.995 (dotted line), 1.98 (dashed 
line), 1.5 (long-dashcd line) and (solid line) for the pure gas collapse with parameters M = 10^^ Mq, 
Zc = 2 and To = 0.014 K . 
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Fig. 2.— Over-density profiles versus tlie self-similar radius A (see Section 3.4) obtained at different redsliifts 
for the pure dark matter case. We used scale invariant initial conditions with e = 0.2, 0.5, 0.8 and 1. The 
self similarity is perfectly recovered, as well as the correct power law. The position of the dark matter "shock 
front" (arrows) is exactly at the radius of the first caustic calculated by Fillmore & Goldreich (1984). 
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Fig. 3. — Over-density profiles versus radius obtained at redshifts z = 1.995 (dotted line), 1.98 (dashed 
line), 1.5 (long-dashed line) and (solid line) for the pure dark matter case. We used a single Fourier 
spherical mode as initial perturbation. Before the dark matter "shock front" (arrows), the pre-coUapse slope 
n = —12/7 is recovered. After that, the corresponding self similarity solution of Fillmore & Goldreich (1984) 
is obtained with a slope close to n = —2. 
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Fig. 4. — Dynamical structure of a spherical cluster at z = with various fls- (a) Gas over-density profiles 
for f^B = 1 (solid line), 0.99 (dotted lino), 0.9 (dashed line) and 0.01 (long-dashed hne). We also plot for 
comparison the dark matter over-density profile obtained in the fis = case (upper solid line), (b) Dark 
matter over-density profile for the same values of Os. Note that all profiles are similar to the Os = case, 
(c) Gas temperature for Qb = 1, 0.99, 0.9 and 0.01. (d) Gas mass fraction versus total mass for fig = 1, 
0.99, 0.9 and 0.01. 
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Fig. 5. — Dynamical structure of a spherical cluster (M = 10^^ Mq and Zc = 2) obtained at z = for various 
initial entropy levels. The parameters are (see text) Tbg = 10^ K (dotted line), 10^ K (dashed line) and 10^ 
K (solid line). 



